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Abstract We propose an experiment to measure the 
slow log(A'^) convergence to mean-field theory (MFT) around 
a dynamical instability. Using a density matrix formalism, we 
derive equations of motion which go beyond MFT and provide 
accurate predictions for the quantum break-time. The leading 
quantum corrections appear as decoherence of the reduced 
single-particle quantum state. 

A Bose-Einstein condensate is described in mean field 
theory (MFT) by a c-number macroscopic wave func- 
tion, obeying the Gross-Pitaevskii non-linear Schrodinger 
equation. MFT is closely analogous to the semiclassi- 
cal approximation of single-particle quantum mechanics, 
with the inverse square root of the number N of particles 
in the condensate playing the role of ?i as a perturbative 
parameter. Since in current experimental condensates N 
is indeed large, it is generally difficult to see qualitatively 
significant quantum corrections to MFT. In the vicinity 
of a dynamical instability in MFT, however, quantum 
corrections appear on timescales that grow only loga- 
rithmically with N . In this paper we propose an exper- 
iment to detect such quantum corrections, and present 
a simple theory to predict them. We show that, as the 
Gross-Pitaevskii classical limit of a condensate resembles 
single-particle quantum mechanics, so the leading quan- 
tum corrections appear in the single-particle picture as 
decoherence. 

We will consider a condensate in which particles can 
only effectively populate two second-quantized modes. 
This model can be realized with a condensate in a double 
well potential , or with an effectively two-component 
spinor condensate whose internal state remains uni- 
form in space. Such uniformity can be ensured, to a 
good approximation, by confining a very cold condensate 
within a size much smaller than [n\^Jaiia22 — ai2|]~^^^j 
where n is the mean total density and is the s-wave 
scattering length between atoms in internal states i and 
j. The kinetic energy of spin non- uniformity then ensures 
that the spatial state of the condensate will adiabatically 
follow its internal state, which will evolve on slower time 
scales. Dynamical instabilities to phase separation are 
also frustrated in this regime, which since for available 
alkali gases all differ only by a few percent, could 
be reached with small condensates [N < 10^) in weak, 
nearly spherical traps {uj < 100 Hz). Stronger or less 
isotropic traps reach the two-mode regime at smaller N. 

In the double well realization, the nonlinear interac- 
tion may be taken to affect only atoms within the same 
well. In this case single-particle tunneling provides a lin- 
ear coupling between the two modes, which can in princi- 



ple be tuned over a wide range of strengths. Two internal 
states may be coupled by a near-resonant radiation field 
If collisions do not change spin states, there is also 
a simple nonlinear interaction in the internal realization. 
In either case the total number operator commutes with 
the Hamiltonian, and may be replaced with the c-number 
N . Discarding c-number terms, we may therefore write 
the two-mode Hamiltonian 
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{a\)ai - {al)a2 , (1) 



where lo is the coupling strength between the two con- 
densate modes, rj is the two-body interaction strength, 
and ai,a|,a2,a| are particle annihilation and creation 
operators for the two modes. We will take rj and w to be 
positive, since the relative phase between the two modes 
may be re-defined arbitrarily, and since without dissipa- 
tion the overall sign of H is insignificant. 

Instead of considering the evolution of a-,- and its ex- 
pectation value in a symmetry-breaking ansatz, we will 
examine the evolution of the directly observable quan- 
tities ajaj, whose expectation values define the reduced 
single particle density matrix (SPDM) Rij = {a\aj)/N. 
It is convenient to introduce the Bloch representation, by 
defining the angular momentum operators, 

p _ a\a2 + 0.201 f _ a\a2 — d\ai 
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^ _ d\ai — d\d2 
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The Hamiltonian Eq. (n|) then assumes the form. 



(2) 



(3) 



and the Heisenberg equations of motion for the three an- 
gular momentum operators of Eq. (0) read 



dt 
d 



\{L.,L, + L,Lx) , (4) 



dt 



The mean-field equations for the SPDM in the two- 
mode model may be obtained, without invoking U(l) 
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symmetry breaking, by approximating second-order ex- 
pectation values (LiLj) as products of the first order ex- 
pectation values (Li) and (Lj): 

(L.4-)«(L,)(4) . (5) 

Defining the single-particle Bloch vector 

s = {S,,Sy,S,) = n = 7JN/2 and 

using Eq. (||), we obtain the nonlinear Bloch equations 

Sr — kSzS^i 



Jy , 



(6) 



Mean field trajectories s(t) at four different k/uj ratios are 
plotted in Fig. 1. The norm of s is conserved in MFT, 
and so for an initially pure SPDM, are equivalent to 
the two-mode Gross-Pitaevskii equation 0] . 
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If we assume that the condensate remains mildly frag- 
mented, so that the two eigenvalues of Rij are / and 
1 — / for small /, we can take Li — Li + 5 Li, where 
the c- number Li is 0{N) and, throughout the Hilbert 
subspace through which the system will evolve, all ma- 
trix elements of 5Li are no greater than 0{Ny/f ). The 
second order moments 



LjLi 



2{U){L, 



(7) 



will then be of order /. We can retain these, and so 
improve on MFT, if we truncate the BBGKY hierarchy of 
expectation value equations of motion at one level deeper: 
we eliminate the approximation (|^), and instead impose 

{LiLjLk) w (LiLj){Lk) + {Li){LjLk) + (LiLk){Lj) 

^2{U){L,){U) . (8) 

This approximation is accurate to within a factor 1 + 
better than (|) by one factor of f^'^. Succes- 
sively deeper truncations of the hierarchy yield system- 
atically better approximations as long as / is small. 

Applying (|^) and (^) to we obtain the following set 
of equations, in which the mean field Bloch vector drives 
the fluctuations Ay, and is in turn subject to backreac- 
tion from them: 




FIG. 1. Mean-field trajectories at {&)k = 0, (b)K = 1.02aj, 
(c)k = 2uj, and (d)K = 2Qu}. 

The nonlinear Bloch equations (^) depict a competi- 
tion between linear Rabi oscillations in the S'yS^-plane 
and nonlinear oscillations in the S'^S'^-plane. For a non- 
interacting condensate (Fig. la) the trajectories on the 
Bloch sphere are circles about the Sx axis, correspond- 
ing to harmonic Rabi oscillations. As k increases the 
oscillations become more anharmonic until, above the 
critical value k = uj (Fig. lb), the stationary point 
s = (—1, 0, 0) becomes dynamically unstable, and macro- 
scopic self-trapping can occur (oscillations with a non- 
vanishing time averaged population imbalance {Sz)t 7^ 0) 
In the vicinity of the dynamically unstable point, 
MFT will break down on a time scale only logarithmic 
in N , and so an improved theory is desirable. 



Sx I^SzSy 2^ 

Sy ~ i^Sz 4" I^SzSx 
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Sz — —ojS^i 



A,. = 

Axy = 
^xx — 

^yy ~ 

Azz = 



/^*S^ A,, 



I^Sy^Z 



■yz 

(^i^zz ~ Ayy) + kSz^xz + 
{UJ + kSx)Axz - KSyAyz + 

lnSyAxz ^KjSzAxy 
2{UJ + kSx)\jz + 2KSzAxy 
-2LUAyz . 



kSxAzz 
kSz{Ax 



(9) 
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In what follows, we will refer to (||) as evolution under 
'Bogoliubov backreaction' (BBR). 

We note that (||) are actually identical to the equations 
of motion one would obtain, for the same quantities, us- 
ing the Hartree-Fock-Bogoliubov Gaussian ansatz. And 
if the second-order moments Ay may initially be factor- 
ized as Ay — SiSj = x,y,z), then the factorization 
persists and the time evolution of Sx,Sy, and Sz is equiv- 
alent to that of perturbations of the mean-field equations 



Sx = -K{SzSy + SySz) , 

Sy = UjSz + K.{SzSx + SxSz) , 

Sr ~ —i^S„ . 



(10) 



2 



Thus our equations for Ay are in a sense equivalent to the 
usual Bogoliubov equations. The quantitative advantage 
of our approach therefore lies entirely in the wider range 
of initial conditions that it admits, which may more accu- 
rately represent the exact initial conditions. For instance, 
a Gaussian approximation will have A^x = in the 

ground state, where in fact A^x = 0{N^^). This leads 
to an error of order N~^^'^ in the Josephson frequency 
computed by linearizing (^) around the ground state, 
which is problematic because one would hope that the 
Gaussian backreaction result would be accurate at this 
order. Our SPDM approach avoids this problem, which is 
presumably the two-mode analogue of the Hartree-Fock- 
Bogoliubov spectral gap . 

For finite motion of the Bloch vector, our formalism 
offers an efficient method to depict the back-reaction of 
the Bogoliubov equations on the mean field equations via 
the couphng terms —KAyz/2 and kAxz/2 in (||). Because 
in general Ayz {t) ^ A^z (t) , this back-reaction has the ef- 
fect of breaking the unitarity of the mean- field dynamics. 
Consequently, the BBR trajectories are no longer con- 
fined to the surface of the Bloch sphere, but penetrate to 
the interior (representing mixed-state Rij , with two non- 
zero eigenvalues). Thus although decoherence is gener- 
ally considered as suppressing quantum effects, decoher- 
ence of the single particle quantum state of a condensate 
is itself the leading quantum correction (due to interpar- 
ticle entanglement) to the effectively classical MFT. 




FIG. 2. Mean field (■■•), Bogoliubov back-reaction ( ) 

and exact 50 particles ( ) trajectories starting with all par- 
ticles in one mode, at k = 2w. 

In order to confirm this decoherence, and demonstrate 
how the BBR equations (^) improve on the mean-field 
equations (|^), we compare the trajectories obtained by 
these two formalisms to exact quantum trajectories, ob- 
tained by fixing the total number of particles = 50 



and solving Eq. (|^) numerically, using a Runge-Kutta 
algorithm. The results for an initial state where all par- 
ticles are in one of the modes (corresponding to the ini- 
tial conditions Sx = Sy = 0,Sz = —l,Axx = —Ayy = 
2/N,Axy = Axz = Ayz = Azz = 0) and k = 2uj are 
shown in Fig. 2. The MFT trajectory passes through the 
dynamically unstable point s = (—1,0,0). Consequently, 
the quantum trajectory sharply breaks away from the 
MFT trajectory as it approaches this point, entering the 
Bloch sphere interior. While still periodic on a much 
shorter time scale than the exact evolution, the BBR evo- 
lution (dashed curve) provides an excellent prediction of 
the time at which the break from MFT takes place (the 
'quantum break time'). 
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FIG. 3. Growth of the von Neumann entropy S of the quan- 
tum reduced single-particle density operator, at k = 2u}, for 

N=10 (- • -), 20 (• • •), 40 ( ), 80 ( ), 160 (- . - . -), 

and 320 ( ) particles. Initial conditions are the same as 

in Fig. 2. 

The quantum break time near a dynamical instability 
is expected to grow logarithmically with N. In Fig. 3 we 
plot the von Neumann entropy 

5 = -li„["^l-T"""7-l'T"""'] (11) 

of the exact reduced single-particle density operator, as a 
function of the rescaled time ujt with N=10, 20,40, 80, 160, 
and 320 particles, for the same initial conditions as in 
Fig. 2. Since the MFT entropy is always zero, S serves 
as a measure of convergence. The quantum break time 
is clearly evident, and indeed increases as log(A^). The 
single-particle entropy is measurable, in the internal state 
realization of our model, by applying a fast Rabi pulse 
and measuring the amplitude of the ensuing Rabi oscil- 
lations, which is proportional to the Bloch vector length 
|s|. (Successive measurements with Rabi rotations about 
different axes, i.e. by two resonant pulses differing by 
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a phase of 7r/2, will control for the dependence on the 
angle of s). In a double well realization, one could deter- 
mine the single-particle entropy by lowering the potential 
barrier, at a moment when the populations on each side 
were predicted to be equal, to let the two parts of the 
condensate interfere. The fringe visibility would then be 
proportional to |s|. Since |s| at a fixed time depends ex- 
ponentially on 77, such experiments could also potentially 
be used to measure scattering lengths. 
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FIG. 4. Time at which S reaches 0.2 as a function of the 
particle number A'', according to the BBR equations mod- 
ified to include thermal phase-diffusion. Four different values 

of r are shown: F = ( ), F = IQ-'^uj (■ ■ ■),V = IO'^lj 

( ), and F — IQ'^uj (— ■ Exact quantum results are 

presented for F = (circles) and F = IO~^tJ (squares). Initial 
conditions, k and uj are the same as in Figs. 2 and 3. 

Decoherence of quantum systems coupled to reservoirs 
shows similar behaviour to Fig. (3). The entropy of a dy- 
namically unstable quantum system coupled to a reser- 
voir jlj], or of a stable system coupled to a dynami- 
cally unstable reservoir, is predicted to grow linearly with 
time, at a rate independent of the system-reservoir cou- 
pling, after an onset time proportional to the logarithm 
of the coupling. This shows that one can really consider 
the Bogoliubov fluctuations as a reservoir jl^, coupled 
to the mean field with a strength proportional to 1/N. 
But one can also consider decoherence due to a genuine 
reservoir (unobserved degrees of freedom, as opposed to 
unobserved higher moments). For example, thermal par- 
ticles scattering off the condensate mean field will cause 
phase diffusion jl^ at a rate F which may be estimated 
in quantum kinetic theory as proportional to the ther- 
mal cloud temperature. For internal states not entan- 
gled with the condensate spatial state, F may be as low 
as 10~^ Hz under the coldest experimental conditions, 
whereas for a double well the rate may reach 10~^ Hz. 
Further sources of decoherence may be described phe- 
nomenologically with a larger F. Evolving the full A^- 
particle density matrix under the appropriate quantum 



kinetic master equation Q , we again solve for s either nu- 
merically or in BBR approximation. In Fig. 4 we show 
the time at which the entropy reaches a given value, as a 
function of the number of particles, for various F, accord- 
ing to the modified BBR equations. The exact quantum 
results (limited by computation power to ~ 10^) are 
presented for the two limiting values of F, showing excel- 
lent agreement with the BBR predictions. These results 
can of course be interpreted as a quantum saturation of 
the dephasing rate at low temperature. 
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